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ABSTRACT 

We have constructed the first all-sky map of the thermal Sunyaev-Zeldovich (tSZ) effect by applying specifically tailored component separation 
algorithms to the 100 to 857 GHz frequency channel maps from the Planck survey. These maps show an obvious galaxy cluster tSZ signal that is 
well matched with blindly detected clusters in the Planck SZ catalogue. To characterize the signal in the tSZ map we have computed its angular 
power spectrum. At large angular scales (I < 60), the major foreground contaminant is the diffuse thermal dust emission. At small angular 
scales (€ > 500) the clustered Cosmic Infrared Background (CIB) and residual point sources are the major contaminants. These foregrounds are 
carefully modelled and subtracted. We measure the tSZ power spectrum in angular scales, 0.17° < 6 < 3.0°, that were previously unexplored. The 
measured tSZ power spectrum is consistent with that expected from the Planck catalogue of SZ sources, with additional clear evidence of signal 
from unresolved clusters and, potentially, diffuse warm baryons. We use the tSZ power spectrum to obtain the following cosmological constraints: 
cr 8 (O m /0.28) 3 2/8 1 = 0.784±0.016 (68 % C.L.). Marginalized band-powers of the Planck tSZ power spectrum and the best-fit model are given. The 
non-Gaussianity of the Compton parameter map is further characterized by computing its ID probability distribution function and its bispectrum. 
These are used to place additional independent constraints on tr 8 . 

Key words, cosmological parameters - large-scale structure of Universe - Galaxies: clusters: general 
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1. Introduction 

This paper, one of a set associated with the 2013 release of data 
from the Planck 1 mission (Planck Collaboration I 2013), de- 
scribes the construction of a Compton parameter map and its 
angular power spectrum and high order statistics. 

The thermal Sunyaev-Zeldovich (tSZ) effect (Sunyaev & 
Zeldovich 1972), produced by the inverse Compton scattering 
of CMB photons by hot electrons along the line of sight, has 
proved to be a major tool to the study of the physics of clusters 
of galaxies as well as structure formation in the Universe. In par- 
ticular, tSZ-selected catalogues of clusters of galaxies have been 
provided by various experiments including the Planck satellite 
(Planck Collaboration VIII 2011; Planck Collaboration XXIX 
2013), the Atacama Cosmology Telescope (ACT, Hasselfield 
et al. 2013) and the South Pole Telescope (SPT, Reichardt 
et al. 2013). From these catalogues, and their associated sky 
surveys, a wealth of studies have been performed both on the 
physics of clusters of galaxies (Planck Collaboration XII 2011; 
Planck Collaboration XI 2011; Planck Collaboration X 2011) 
and on their cosmological interpretation (Planck Collaboration 
XX 2013; Benson et al. 2013; Das et al. 2013; Wilson et al. 2012; 
Mak & Pierpaoli 2012). 

The study of number counts and their evolution with red- 
shift using tSZ detected clusters of galaxies has been recog- 
nized as an important cosmological test (Carlstrom et al. 2002; 
Dunkley et al. 2013; Benson et al. 2013; Planck Collaboration 
XX 2013). As a complement, the measurement of the tSZ ef- 
fect power spectrum has been proposed by Komatsu & Seljak 
(2002). One advantage of using the tSZ angular power spectrum 
over cluster counts is that no explicit measurement of cluster 
masses is required. However, drawbacks of using the tSZ angu- 
lar power spectrum include potential contamination from point 
sources (Rubino-Martm & Sunyaev 2003; Taburet et al. 2010) 
and other foregrounds. Also, lower mass, and therefore fainter, 
clusters, that may not be significantly detected as individual ob- 
jects contribute to this statistical signal (Battaglia et al. 2010; 
Shaw et al. 2010). 

To date, indirect measurements of the tSZ power spectrum 
are only available from high resolution CMB oriented experi- 
ments like ACT (Sievers et al. 2013) and SPT (Reichardt et al. 
2012). In these studies, constraints on the amplitude of tSZ 
power spectrum at ( — 3000 are obtained by fitting a tSZ tem- 
plate in addition to other components (i.e., CMB, radio and in- 
frared point-source and clustered Cosmic Infrared Background 
(CIB)) to the measured total power spectrum. These constraints 
are obtained at angular scales where the tSZ signal dominates 
over the CMB, but at these same scales the contamination from 
point sources and clustered CIB is important and may affect the 
measured tSZ signal. Moreover, the scales probed are particu- 
larly sensitive to the uncertainties in modelling the intra-cluster 
medium (ICM) over a broad range of masses and redshifts, 
and at large cluster-centric radii (Battaglia et al. 2010). Recent 
work, using hydrodynamical simulations (Battaglia et al. 2010; 
Battaglia et al. 2012) N-body simulations plus semi-analytic gas 
models (Trac et al. 201 1) and analytic models (Shaw et al. 2010), 
have significantly reduced the tension between the observed 
and predicted values. However the distribution of amplitudes 



1 Planck (http://www.esa.int/Planck) is a project of the 
European Space Agency (ESA) with instruments provided by two sci- 
entific consortia funded by ESA member states (in particular the lead 
countries France and Italy), with contributions from NASA (USA) and 
telescope reflectors provided by a collaboration between ESA and a sci- 
entific consortium led and funded by Denmark. 



Table 1. tSZ Compton parameter y conversion factors to CMB 
temperature units and the FWHM of the beam of the Planck 
channel maps. 



Frequency 


7CMB g(v) 


FWHM 


[GHz] 


[Kcmb] 


['] 


100 


-4.031 


9.66 


143 


-2.785 


7.27 


217 


0.187 


5.01 


353 


6.205 


4.86 


545 


14.455 


4.84 


857 


26.335 


4.63 



between different models and simulations is still significantly 
larger than the measurement errors, degrading the constraints 
that can be placed on cosmological parameters with these meth- 
ods (Dunkley et al. 2013; Reichardt et al. 2013). 

In addition to the power spectrum, and as pointed out in 
(Rubino-Martm & Sunyaev 2003), the skewness (or equiva- 
lently, the bispectrum) of the tSZ signal is a powerful and inde- 
pendent tool to study and to isolate the signal of clusters, sepa- 
rating it from the contribution of radio and IR sources. Recently, 
Bhattacharya et al. (2012) showed that the bispectrum of the 
tSZ effect signal is dominated by massive clusters at interme- 
diate redshift for which high-precision X-ray observations ex- 
ist. This contrasts with the power spectrum, where the signal 
mainly comes from the lower mass and higher redshift groups 
and clusters (e.g., Trac et al. 2011). The theoretical uncertainty 
in the tSZ bispectrum is thus expected to be significantly smaller 
than that of the SZ power spectrum. Combined measurements of 
the power spectrum and the bispectrum can thus be used to dis- 
tinguish the contribution to the power spectrum from different 
cluster masses and redshift ranges. The bispectrum amplitude 
scales as o~ l & °~ 12 (Bhattacharya et al. 2012). Similarly, Wilson 
et al. (2012) used the unnormalized skewness of the tSZ fluctu- 
ations, ^r 3 (n)^, which scales as cr l & 01 ~ 11 1 , to obtain an indepen- 
dent determination of <x 8 . 

Thanks to its all sky coverage and unprecedented wide fre- 
quency range, Planck has the unique capability to produce an 
all-sky Compton parameter (y) map and an accurate measure- 
ment of the tSZ power spectrum at intermediate and large an- 
gular scales, for which the tSZ fluctuations are almost insensi- 
tive to the cluster core physics. The Planck Compton parame- 
ter map also offers the possibility to study the properties of the 
non-Gaussianity of the tSZ signal using higher order statistical 
estimators as the skewness and the bispectrum. In this paper we 
derive such a tSZ all-sky map from the individual Planck fre- 
quency maps and compute its power spectrum, its ID probability 
density function (ID PDF), and the associated bispectrum. 

The paper is structured as follows. Sect. 2 describes the 
Planck data used to compute the tSZ all-sky map and the simula- 
tions used to characterize it. We discuss details of the modelling 
of the tSZ effect power spectrum and bispectrum in Sect. 3. In 
Sect. 4 we present the Planck all-sky Compton parameter map. 
Sect. 5 describes the power spectrum analysis. Cross-checks us- 
ing high-order statistics are presented in Sect. 6. Cosmological 
interpretation of the results is discussed in Sect. 7 and we present 
our conclusions in Sect. 8. 
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2. Data and simulations 

2.1. The Planck Data 

This paper is based on the first 15.5 month survey mission 
corresponding to two full-sky surveys. We refer to Planck 
Collaboration II (2013) and Planck Collaboration VI (2013) for 
the generic scheme of Time Ordered Information (TOI) process- 
ing and map-making, as well as for the technical characteris- 
tics of the Planck frequency maps. The Planck channel maps 
are provided in HEALPix (Gorski et al. 2005) N s i de = 2048 
at full resolution. An error map is associated with each chan- 
nel map and is obtained from the difference of the first half and 
second half of the survey rings for a given pointing position of 
the satellite spin axis. The resulting maps, called null maps, are 
mainly free from astrophysical emission and they are a good rep- 
resentation of the statistical instrumental noise. Null maps have 
also been used to estimate the noise in the final Compton pa- 
rameter maps. Here we approximate the Planck beams by effec- 
tive circular Gaussians (Planck Collaboration IV 2013; Planck 
Collaboration VII 2013) The FWHM for each frequency chan- 
nel are given in Table 1. Although tests have been performed 
using both LFI and HFI channel maps, the work presented here 
will refer mainly to results using the HFI data only. 

2.2. ffp6 Simulations 

We also use simulated Planck frequency maps obtained from 
the so-called Full Focal Plane (FFP6) simulations (described 
in Planck Collaboration ES 2013)). These simulations include 
the most relevant sky components at microwave and millime- 
ter frequencies: CMB; thermal SZ effect; diffuse Galactic emis- 
sions (synchrotron, free-free, thermal and spinning dust); radio 
and infrared point sources, and the clustered Cosmic Infrared 
Background (CIB) based on foreground models from the Planck 
Sky Model Delabrouille et al. (PSM, 2012). The simulated tSZ 
signal was constructed using hydrodynamical simulations of 
clusters of galaxies up to redshift 0.3, completed with pressure 
profile-based simulations of individual clusters of galaxies ran- 
domly drawn on the sky. The noise in the maps was obtained 
from Gaussian random realizations of noise in the time domain 
and therefore accounts for noise inhomogeneities in the maps. 



3. Modelling the tSZ effect 

The thermal SZ Compton parameter in a given direction, n, is 



v( 



n) = J n e 



-CTj ds 



(1) 



where ds is the distance along the line-of-sight, n, and « e and T e 
are the electron number density and temperature, respectively. 

In units of CMB temperature the contribution of the tSZ ef- 
fect to the Planck maps for a given observation frequency v is 



Ar 



! CMB 



; g(y) y- 



(2) 



Neglecting relativistic corrections we have g(v) = 
coth(|) - 4) , with x = hv/(k B TcMB)- Table 1 shows 
the Compton parameter to CMB temperature, Kcmb> conversion 
factors for each frequency channel after integrating within the 
bandwidth. 



3.1. tSZ power spectrum 

Decomposing the map in spherical harmonics we obtain 



y(n) = 2 

(in 



ytm Y{ m (n), 



(3) 



Thus, the angular power spectrum of the Compton parameter 
map is 



C 



tSZ 



1 y * 
2 / + 1 2-i ytmy ^ 



(4) 



Note that Cf z is an a-dimensional quantity as y. 

To model the tSZ power spectrum we consider a two halo 
model to account for intra-halo and inter-halo correlations 



£-iSZ £lhalo _j_ £2halos 



(5) 



The 1 halo term, also known as the Poissonian contribution, can 
be computed by summing the square of the Fourier transform 
of the projected SZ profile, weighted by the number density of 
clusters of a given mass and redshift (Komatsu & Seljak 2002): 



Jo dzdQ, Jm 



dM d -^ mM , z) \\ (6) 
dM 



where dV c /(dzdQ) is the comoving volume per unit redshift and 
solid angle and n(M, z)dM dV c /(dzdQ) is the probability of hav- 
ing a galaxy cluster of mass M at a redshift z in the direction 
dQ. The quantity JV = ye (M, z) is the 2D Fourier transform on 
the sphere of the 3D radial profile of the Compton ^-parameter 
of individual clusters, 



l{ \ m e c z 



)I 



dxJP e {M,z,x)^J^ (7 ) 



where x = r/r s , £ s = DpXz)lr s , r s is the scale radius of the 
3D pressure profile, D^(z) is the angular diameter distance to 
redshift z and P e is the electron pressure profile. 

The two-halo term is obtained by computing the correlation be- 
tween two different halos (Komatsu & Kitayama 1999; Diego & 
Majumdar 2004; Taburet et al. 2011) 



^2halos 



-r 

(I 



dz 



dV c 
dzdQ, 



dn(M, z) 
dM \ ' \y c (M,z)\ B(M,z) 
dM 



P(k,z), (8) 



where P(k,z) is the 3D matter power spectrum at redshift z. 
Here B(M,z) is the time-dependent linear bias factor that re- 
lates the matter power spectrum, P(k,z), to the power spec- 
trum of the cluster correlation function. Following Komatsu & 
Kitayama (1999, see also Mo & White 1996) we adopt B(M, z) = 
1 +(v\M,z)- l)/<5 c (z), where v(M,z) = 6 c (M)/D(z)cr(M), cr(M) 
is the present-day rms mass fluctuation, D(z) is the linear growth 
factor, and 6 c (z) is the threshold over-density of spherical col- 
lapse. 

Finally, we compute the tSZ power spectrum using the 
Tinker et al. (2008) mass function dn{M,z)/dM including an 
observed-to-true mass bias of 20%, as discussed in detail in 
Planck Collaboration XX (2013), and we model the SZ Compton 
parameter using the pressure profile of Arnaud et al. (2010). This 
approach is consistent with the ingredients of the cluster number 
count analysis in Planck Collaboration XX (2013). 
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3.2. Nth moment of the tSZ field 



4.1 .1 . MILCA 



To calculate the Mh moment of the tSZ field, we assume to 
first order that the distribution of clusters on the sky can be ade- 
quately described by a Poisson distribution corresponding to the 
one halo term. We neglect the contribution due to clustering be- 
tween clusters and their overlap (Komatsu & Kitayama 1999). 
The Mh moment is then given by (Wilson et al. 2012) 



f W dz dVc f 
Jo dzdD, J Mn 



dM 



dn(M, z) 
dM 



J d 2 0y(6,M,zf, 



(9) 



where y(6, M, z) is the integrated Compton parameter along the 
line-of-sight for a cluster of mass M at redshift z. 

3.3. Bispectrum 

The angular bispectrum, analogous to the 3-point correlation 
function in harmonic space, is the lowest order indicator of the 
non-Gaussianity of a field. It is given by 



B 



= (ye i mje 1 m 2 ye i m i ) , 



(10) 



where the angle-averaged quantity in the full-sky limit can be 
written as 



*(W 3 )= X (£ 

m. m-i m-5 * 



ni2 mi 



nm l m 2 m } 
W3 ' 



(ID 



which has to satisfy the conditions: m.\ +mi+m 3 = 0; l\ +^2+^3 
even : and \£j - t\ < 4 < (i + (j for the Wigner 3 j function in 
brackets. For illustration we compute the bispectrum assuming a 
Poissonian distribution, given by (Bhattacharya et al. 2012): 



dV c 







Jo dzdQ. J M . 



4n 



dM dH{ ^ M Z) y ti (M, z)MM, Z)%(M, zJ12) 



4. The reconstructed all-sky tSZ map 

4. 1. Reconstruction methods 

The contribution of the tSZ effect in the Planck frequency maps 
is subdominant with respect to the CMB and other foreground 
emissions. Furthermore, the tSZ effect from galaxy clusters is 
spatially localized and leads to a highly non-Gaussian signal 
with respect to that from the CMB. CMB-oriented component- 
separation methods (Planck Collaboration XII 2013) are not 
optimized to recover the tSZ signal. We therefore need to 
use specifically tailored component separation algorithms that 
are able to reconstruct the tSZ signal from the Planck fre- 
quency channel maps. These optimized all-sky component sep- 
aration techniques rely on the spatial localization of the dif- 
ferent astrophysical components and on their spectral diversity 
to separate them. We present,in the following, the results of 
two algorithms, MILCA (Modified Internal Linear Combination 
Algorithm, Hurier et al. 2010) and NILC (Needlet Iindependent 
Linear Combination, Remazeilles et al. 201 1). Both are based on 
the well known Internal Linear Combination (ILC) approach that 
searches for the linear combination of the input maps that mini- 
mizes the variance of the final reconstructed map under the con- 
straint of offering unit gain to the component of interest (here the 
tSZ effect, whose frequency dependence is known). Both have 
been extensively tested on simulated Planck data. 



MILCA (Hurier et al. 2010) uses two constraints: preservation of 
the tSZ signal, assuming the tSZ spectral signature; and removal 
of the CMB contamination in the final SZ map, making use of 
the well known spectrum of the CMB. In addition, to compute 
the weights of the linear combination, we have used the extra 
degrees of freedom in the linear system to minimize residuals 
from other components (two degrees) and from the noise (two 
degrees). The noise covariance matrix was estimated from the 
null maps described in Section 2. 1 . To improve the efficiency of 
the MILCA algorithm, weights are allowed to vary as a function 
of multipole £ , and are computed independently on different sky 
regions. We have used 11 filters in ( space with an overall trans- 
mission of one, except for t < 8. For these large angular scales 
we have used a Gaussian filter to reduce foreground contami- 
nation. The size of the independent sky regions was adapted to 
the multipole range to ensure sufficient spatial localization at the 
required resolution. We used a minimum of 12 regions at low 
resolution and a maximum of 3072 regions at high resolution. 

4.1 .2. NILC 

In the multi-component extensions of NILC Delabrouille et al. 
(2009); Remazeilles et al. (2011), initially developed to extract 
the CMB, the weights for component separation (i.e., covari- 
ances) are computed independently in domains of a needlet de- 
composition (in the spherical wavelet frame). The needlet de- 
composition provides localization of the ILC filters both in pixel 
and in multipole space, allowing us to deal with local contam- 
ination conditions varying both in space and in scale. We im- 
posed constraints to remove the CMB contamination and pre- 
serve the tSZ effect. To avoid strong foreground contamination, 
the Galactic plane was masked before applying NILC to the 
Planck frequency maps. 

In both methods, we mask the brightest regions in the Planck 
857 GHz channel map, corresponding to about 33 % of the sky. 
We use the HFI channel maps from 100 to 857 GHz that are 
convolved to a common resolution of 10'. The 857GHz channel 
map is mainly exploited in the internal linear combination as a 
template to remove the thermal dust emission on large angular 
scales. However it induces significant CIB residuals in the tSZ 
map on small scales. To avoid this contamination, while enabling 
efficient removal of the diffuse thermal dust emission at large 
angular scales, we use the 857GHz channel only for I < 300. 

4.2. Reconstructed Compton parameter y map 

Figure 1 shows the reconstructed Planck all-sky Compton pa- 
rameter map for NILC (top panel) and MILCA (bottom panel). 
For display purposes, the maps are filtered using the procedure 
described in Sect. 6. Clusters appear as positive sources: the 
Coma cluster and Virgo supercluster are clearly visible near the 
north Galactic pole. As mentioned above, the Galactic plane 
is masked in both maps, leaving 67% of sky. Other weaker 
and more compact clusters are visible in the zoomed region of 
the Southern cap shown in the bottom panel of Fig. 2. Strong 
Galactic and extragalactic radio sources show up as negative 
bright spots on the maps and were masked prior to any scientific 
analysis, as discussed below in Sect. 4.3. Residual Galactic con- 
tamination is also visible around the edges of the masked area; 
extra masking was performed to avoid this highly contaminated 
area. The difference of contrast observed between the NILC and 
MILCA maps comes both from differences in the noise and in- 
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NILC tSZ map 





Fig. 1. Reconstructed Planck all-sky Compton parameter maps for NILC (top) and MILCA (bottom) in orthographic projections. 
The difference of contrast observed between the NILC and MILCA maps comes both from differences in the noise and instrumental 
systematic contribution and from the differences in the filtering applied for display purpose to the original Compton parameter 
maps. 
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NILC tSZ map 



MILCA tSZ map 




3.5 ^^^^m 5.0 y x 10° 

(0.0. -45.0) Galactic 



Fig. 2. Zooms into the reconstructed Planck all-sky Compton parameter maps for NILC (left) and MILCA (right) at intermediate 
Galactic latitudes in the southern sky. 



strumental systematic contribution (the NILC map is slightly 
noisier but less affected by systematic effect and foreground 
emissions than the MILCA map, as discussed in Sect. 5.2) and 
from the differences in the filtering applied for display pur- 
pose to the original Compton parameter maps, as discussed in 
Sect. 6.1. 

In addition to the full Compton parameter maps, we also pro- 
duce the so-called FIRST and LAST Compton parameter maps 
from the first and second halves of the survey rings (i.e., pointing 
periods). These maps are used for the power spectrum analysis 
in Sect. 5. 



4.3. Point source contamination and masking 

Point source contamination is an important issue for the cos- 
mological interpretation of the Planck Compton parameter map. 
Radio sources will show up in the reconstructed tSZ maps as 
negative peaks, while infrared sources will show up as posi- 
tive peaks, mimicking the cluster signal. To avoid contamination 
from these sources we introduce a point source mask (PSMASK, 
hereafter). This mask is the union of the individual frequency 
point-source masks discussed in Planck Collaboration XXVIII 
(2013). To test the reliability of this mask we have performed a 
search for negative sources in the Compton parameter maps us- 
ing the MHW2 algorithm (Lopez-Caniego et al. 2006). We have 
found that all resolved radio sources in the Compton parameter 
maps are masked by the PSMASK. For infrared sources, esti- 
mating the efficiency of the masking is hampered by the tSZ 
signal itself. The residual contamination from point sources is 
discussed in Sect. 5.2 and Sect. 6. It is also important to note that 
the PSMASK may also exclude some clusters of galaxies. This is 
particularly true in the case of clusters with strong central radio 
sources, such as the Perseus cluster (see Planck Collaboration 
XXIX 2013, for detailed discussion). 



4.4. tSZ signal from resolved sources 

As a very first validation step of the Compton parameter maps 
we perform a blind search of the SZ signal from resolved sources 
and compare it to the Planck catalogue of SZ sources (Planck 
Collaboration XXIX 2013). The latter comprises 861 confirmed 
clusters out of 1227 cluster candidates and 54 CLASS 1 highly 
reliable candidate clusters. 



4.4.1. Yields 

Two lists of SZ sources above a signal-to-noise threshold of 4.5 
are constructed from both MILCA and NILC all-sky Compton 
parameter maps outside a 33% Galactic mask. The point source 
detections are undertaken using two methods: 

- SMATCH, in which sources are detected using the 
SEXtractor algorithm (Berlin & Arnouts 1996) over the 
whole sky divided into 504 patches. A single frequency 
matched filter (Melin et al. 2006a) is then applied to mea- 
sure the SZ flux density and signal-to-noise using the Arnaud 
et al. (2010) pressure profile. Using this method, we detect 
843 and 872 sources in MILCA and NILC, respectively. 

- MHWS, in which SZ sources are detected in the maps us- 
ing IFCAMEX (MHW2, Gonzalez -Nuevo et al. 2006; Lopez- 
Caniego et al. 2006). The flux density and signal-to-noise are 
then estimated using SEXtractor on 3.65°x3.65° patches. 
We detect 1036 and 1740 sources in MILCA and NILC, re- 
spectively, with this method. 

The difference between the yields of the two methods is un- 
derstandable, as SMATCH is by construction dedicated to the 
search for SZ sources and the precise measurement of their flux 
(including assumptions on the spatial distribution of the SZ sig- 
nal), whereas MHWS targets all types of compact source (includ- 
ing IR and radio sources) and uses a more "generic" flux estima- 
tion procedure. 
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Fig. 3. Comparison of the measurement of Fjrsoo- Left: between the values derived from the detection methods used to build the 
Planck catalogue of clusters (T 5R5 oo,pszX an d those from the all-sly reconstructed MILCA tSZ map (T 5S5 oo,y-MAp)- Right: between 
the MILCA (T5r5oo,milca) and NILC (Ts«5oo,nilc) all-sky tSZ effect maps. The equality relationship is marked as a dashed black 
line. 



We have compared these two lists of sources with 790 con- 
firmed clusters and CLASS 1 high reliability candidates from 
the Planck catalogue of SZ sources that fall outside the 33% 
Galactic mask. The association is performed on the basis of the 
source positions within a search radius of 10' (the resolution of 
the SZ all-sky maps). We found 583 and 529 matches in the 
MILCA source list with the SMATCH and MHWS methods, re- 
spectively (614 and 414 from the NILC source list). This match 
of about 52 to 77% per cent, respectively, is satisfactory. Indeed, 
as shown in Melin et al. (2012), indirect detection methods based 
on reconstructed >'-maps are less efficient at extracting clusters 
of galaxies than dedicated direct methods such as those used to 
build the Planck catalogue of SZ sources (i.e., MMF1, MMF3 and 
PwS Herranz et al. 2002; Melin et al. 2006b; Carvalho et al. 
201 1; Planck Collaboration XXIX 2013). 

4.4.2. Photometry 

Of more importance than a comparison of yields is the compar- 
ison in terms of photometry. For all-sky map detections that are 
associated with clusters in the Planck SZ catalogue, the SZ flux 
measurement from the all-sky maps correlates very well with the 
maximum likelihood value of the integrated Compton parame- 
ter, T5r 500 2 , provided by the dedicated SZ-detection methods in 
the Planck SZ catalogue. As shown in the top panel of Fig. 3, 
the correlation is very tight, with little dispersion (0. 1 dex) . We 
note that the few points at high T5« 5O0 that lie significantly above 
the one-to-one line are expected; they correspond to nearby and 
extended clusters. On the one hand, the significance of SZ flux 
measurement increases with the flux. On the other hand, the cat- 
alogue detection methods are not optimized for the extraction of 
such extended sources (see Planck Collaboration XXIX 2013, 
for details), therefore they tend to miss part of the SZ flux, which 



2 W500 refers to the radius inside which the mean density is 500 times 
the critical density at the cluster redshift. 



is recovered together with a better estimate of the cluster size 
from the Compton parameter map directly. 

As a sanity check, we have also matched the list of sources 
detected by a given method on both MILCA and NILC maps in 
order to compare the SZ photometry. The bottom panel of Fig. 3 
shows very good agreement between the methods. There is only 
0.07 and 0.01 dex dispersion between them for the SMATCH and 
MHWS extraction methods, respectively. 

Together, these results indicate that we can be confident in 
the fidelity with which the tSZ signal is reconstructed over the 
whole sky by the MILCA and NILC methods. 

5. Angular power spectrum of the reconstructed 
y-map 

5.1. Methodology 

To estimate the power spectrum of the tSZ signal we use the 
XSPECT method (Tristram et al. 2005) initially developed for 
the cross-correlation of independent detector maps. XSPECT 
uses standard MASTER-like techniques (Hivon et al. 2002) to 
correct for the beam convolution and the pixelization, as well 
as the mode-coupling induced by masking foreground contami- 
nated sky regions. 

We apply XSPECT to the FIRST and LASTy-maps obtained 
using NILC and MILCA. We consider the following map pairs: 
the MILCA FIRST and LAST (MILCA F/L); the NILC FIRST 
and LAST (NILC F/L); and the NILC FIRST and MILCA LAST 
(NILC-MILCA F/L), or equivalently the MILCA FIRST and 
NILC LAST (MILCA-NILC F/L). As the noise is decorrelated 
between the map pairs the resulting power spectrum is not biased 
and we preserve the variance. 

In the following, all the spectra will use a common mul- 
tipole binning scheme that was defined in order to minimize 
the correlation between adjacent bins at low multipoles and to 
increase the signal-to-noise at high multipole values. Error bars 
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Fig. 4. Angular power spectrum of the main foreground contri- 
butions as estimated using the FFP 6 simulations. We plot the 
diffuse Galactic emission (blue), clustered CIB (green) and point 
source (cyan) contributions, as well as the tSZ signal (red). The 
solid and dotted lines correspond the NI LC F/L and to the NILC- 
MI LCA F/L cross-power spectra, respectively. For illustration we 
also show the Planck instrumental noise power spectrum (dashed 
black line) in the MI LCA Compton parameter map. 
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Fig. 5. Angular cross-power spectra of the Planck NILC F/L 
reconstructed Compton parameter maps for different Galactic 
masks, removing 30 (cyan), 40 (black points and error bars), 50 
(red), 60 (green) and 70 (blue) % of the sky. 



in the spectrum are computed analytically from the auto-power 
and cross-power spectra of the pairs of maps, as described in 
Tristram et al. (2005). All of our Compton parameter maps 
assume a circular Gaussian beam of 10' FWHM. The additional 
filtering at large angular scales in the MI LCA Compton parame- 
ter maps is also accounted for and deconvolved. 



5.2. Foreground contamination 

The challenge in computing the tSZ power spectrum is to esti- 
mate and minimize foreground contamination. We do not intend 
here to provide a detailed foreground analysis, but rather to iden- 
tify the main foreground contaminants at different multipoles. 
We first identify the dominant foregrounds in the reconstructed 
Compton parameter maps. To do so, we apply to the FFP 6 simu- 
lated maps the linear combination weights of NILC and MILCA 
derived from the real data. In this way we have constructed maps 



of the expected foreground contamination in the final Compton 
parameter maps. 

Figure 4 shows the angular power spectra for these recon- 
structed foreground contamination maps. The PS MASK and a 
conservative common Galactic mask that leaves 50% of the sky 
are used. The Galactic mask is constructed by masking the 50% 
brightest regions of the sky in the 857 GHz intensity map, as de- 
tailed below in Sect. 5.2.1. We show the diffuse Galactic contam- 
ination (blue), the clustered Cosmic Infrared Background con- 
tamination (green), and point source contamination (cyan). We 
consider here the foreground contamination in the cross-power 
spectra of the NILC F/L (dotted lines) and NILC-MILCA F/L 
maps (solid lines). The tSZ power spectrum for the FFP 6 simu- 
lations is plotted in red. For illustration we also show the Planck 
instrumental noise power spectrum (dashed black line) in the 
MILCA Compton parameter map. We clearly observe that, as 
expected, the diffuse Galactic emission, mainly thermal dust, 
dominates the foreground contribution at low multipoles. For 
large multipoles the clustered CIB and point source contribu- 
tions dominate the power spectrum. It is important to notice that 
the tSZ signal dominates the angular power spectrum in the mul- 
tiple range 100 < t < 800. We also note that foreground con- 
tamination differs depending on the reconstruction method. We 
observe that MILCA is more affected by foreground contamina- 
tion. However, we find that at large angular scales the diffuse 
Galactic dust contamination is significantly lower in the NILC- 
MILCA F/L cross-power spectrum than in the NILC F/L cross- 
power spectrum. This indicates that the residual dust contamina- 
tion is not 100 % correlated between the reconstructed MILCA 
and NILC Compton parameter maps. In contrast, the clustered 
CIB and point source contamination levels are similar for the 
two cross-power spectra at high multipoles, indicating that the 
residual contamination is 100% correlated between the MILCA 
and NILC maps. 



5.2.1. Low-multipole contribution 

The diffuse Galactic foreground contribution can be signifi- 
cantly reduced by choosing a more aggressive Galactic mask. 
Assuming that at large angular scales the Compton parameter 
maps are mainly affected by diffuse Galactic dust emission, we 
have tested several Galactic masks by imposing flux cuts on the 
Planck 857 GHz channel intensity map. In particular we con- 
sidered masking out 30, 40, 50, 60 and 70% of the sky. The 
edges of these masks have been apodized to limit ringing effects 
on the reconstruction of the angular power spectrum. Figure 5 
presents the angular cross-power spectrum of the reconstructed 
NILC F/L Compton parameter maps for some of these Galactic 
masks: 30 (cyan); 40 (black); 50 (red); 60 (green); 70 (blue)%; 
and the PSMASK. We find that when masking 40% or more of 
the sky the tSZ angular power spectrum does not change sig- 
nificantly. That is why, conservatively, we select the 50% mask 
(GALMASK50 hereafter) that will be used in the remainder of 
this analysis. 

We checked if the foreground contribution in the recon- 
structed Planck Compton parameter maps also depends on the 
reconstruction method. From the analysis of the FFP 6 simula- 
tions we have found that the contribution from foregrounds in 
the NILC and MILCA Compton parameter maps is not the same, 
and it is not fully correlated. Similar results are found for the 
Planck data. Figure 6 shows the cross-power spectra between 
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Fig. 6. Angular cross-power spectra between the reconstructed 
Planck MILCA F/L (black), the NILC F/L (red) and the NILC- 
MILCA F/L (blue) maps. 
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Fig. 7. NILC F/L cross power spectrum before (black points) 
and after (red points) foreground correction, compared to the 
power spectra of the physically motivated foreground mod- 
els. Specifically we show: clustered CIB (green line); infrared 
sources (cyan line); and radio sources (blue line). The statistical 
(thick line) and total (statistical plus foreground, thin line), un- 
certainties are also shown. We also show the best-fit tSZ power 
spectrum model presented in Sect.7.1 as a solid red line. 

the MILCA F/L maps (black) 3 , the NILC F/L maps (red) and the 
NILC- MILCA F/L maps (blue), as a function of £. We observe 
that the MILCA F/L cross-power spectrum presents a larger am- 
plitude than the NILC F/L cross-power spectrum. This is most 
probably due to a larger foreground contamination in the MILCA 
Compton parameter map. 

In addition, we find that the NILC-MILCA F/L 4 cross-power 
spectrum shows the lowest amplitude at low multipoles ({ < 
100). This is due to a reduction of the dust contamination in the 
cross-correlation of the NILC and MILCA Compton parameter 
maps with respect to the dust contamination in the original maps. 
We also find that the NILC-MILCA F/L lies between the MILCA 
F/L and NILC F-L cross-power spectra at high multipoles. This 
can be explained by the differences in the clustered CIB contam- 
ination in the MILCA and NILC Compton parameter maps. An 



3 The excess of power at low E observed in the MILCA F/L maps 
angular cross-power spectrum is due to the deconvolution from the extra 
low-multipole filtering in the MILCA maps, discussed in Sect. 4.1 

4 And equivalently MILCA-NILC F/L that is not shown in the Figure. 



accurate model of the clustered CIB power spectrum is available. 
However, this is not the case for the dust contamination power 
spectrum, thus we restrict the power spectrum analysis presented 
in Sect. 7.1 to t> 60. 

Hereafter, we will consider the NILC F/L cross-power spec- 
trum as a baseline for cosmological analysis and the NILC- 
MILCA F/L cross-power spectrum will be used to cross-check 
the results. 



5.2.2. High-multipole contribution 

The high-/ 1 contamination from clustered CIB and point sources 
affects the measurement of the tSZ spectrum and its cosmolog- 
ical interpretation. Realistic models fitted to the Planck data are 
thus needed. We take advantage of the capability of the Planck 
survey to measure and constrain these foreground emissions and 
use the outputs of Planck Collaboration XVIII (201 1) and Planck 
Collaboration (2013) for the clustered CIB modelling. For the 
six Planck HFI frequencies considered in this paper, the clus- 
tered CIB model consists of six auto-power spectra and 24 cross- 
power spectra. For frequencies above 217 GHz, these spectra 
are fitted in Planck Collaboration (2013) to the measured CIB 
consistently with Planck Collaboration XVIII (201 1). The model 
is extrapolated at 100 and 143 GHz following Bethermin et al. 
(2012) and Planck Collaboration XVIII (2011). The uncertain- 
ties in the clustered-CIB model are mainly due to the cross- 
correlation coefficients that relate the cross-power spectra to the 
auto-power spectra. Following Planck Collaboration (2013) we 
consider 5% global uncertainties on those coefficients. 

We use the Bethermin et al. (2012) model to compute the 
star-forming dusty galaxy contribution. Finally, we use the the 
Tucci et al. (2011) model, fitted to the Planck ERCSC (Planck 
Collaboration Int. VII 2013), for extragalactic radio sources. 
Notice that these models are also used for the study of the clus- 
tered CIB with Planck (Planck Collaboration 2013). 

We now estimate the residual power spectrum in the j-map 
after component separation. We apply the MILCA or NILC 
weights to Gaussian-realization maps drawn using the cross- 
and auto-spectra of each component at the six Planck HFI 
frequencies. The residual power spectrum in the y-map can be 
also estimated in the spherical harmonic domain, as detailed in 
Appendix A. We have tested the consistency between the two 
approaches and we give here results for a map-based estimate 
using a total of 50 all-sky simulations for each of the foreground 
components. Specific simulations, varying the foreground 
models, were also performed to propagate the 5% global un- 
certainties of the model-coefficients into the estimated residual 
power spectrum. We find a 50% uncertainty in the amplitude 
of each residual spectrum (clustered CIB, star-forming dusty 
galaxies and radio sources) in the y-map. 

Figure 7 shows the NILC F/L cross-power spectrum before 
(black points) and after (red points) foreground correction, us- 
ing the refined foreground models presented above. We also 
show the clustered CIB (green), infrared source (cyan) and ra- 
dio source (blue) power spectrum contributions. 



5.3. Contribution of resolved clusters to the tSZ power 
spectrum 

We simulate the expected Compton parameter map for the de- 
tected and confirmed clusters of galaxies in the Planck cata- 
logue (Planck Collaboration XXIX 2013) from their measured 
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Fig. 8. Comparison of the measured tSZ angular power spectrum 
using the cross of the NILC F/L maps (black) with the expected 
angular power spectrum of the confirmed clusters in the Planck 
Cluster Sample (orange line). In red we plot the NILC F/L cross- 
power-spectrum after masking these clusters. The green points 
correspond to the difference of these two cross-power spectra. 
The cross-power spectrum between the NILC Compton parame- 
ter map and the simulated detected cluster map is shown in blue. 
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Fig. 9. ID PDF for the FFP 6 simulation maps considering the 
MILCA linear combination weights obtained for the real data. 
The tSZ effect (red), diffuse Galactic emission (cyan), clustered 
CIB (blue) and radio source (black) contributions to the ID PDF 
are shown. 



integrated Compton parameter, Fsrsoo- The orange solid line in 
Fig. 8 shows the power spectrum of this simulated map. Figure 8 
also shows the cross-power spectrum of the NILC F/L maps 
(black). In red we plot the cross-power spectrum of the NILC 
F/L maps after masking the confirmed clusters from the PSZ 
catalogue. The green curve corresponds to the difference of the 
two cross power-spectra, with and without masking the clusters. 
It is in good agreement with the modelled power spectrum of 
the confirmed clusters of galaxies. We also compute the cross- 
power spectrum of the simulated cluster map and the Planck re- 
constructed Compton parameter NILC map. This is shown in 
blue in the figure. Here again, the signal is consistent with the 
expected power spectrum of the confirmed Planck clusters of 
galaxies. 

These results show that a significant fraction of the signal in 
the reconstructed Planck Compton parameter maps is due to the 
tSZ effect of detected and confirmed clusters of galaxies, verify- 
ing the SZ nature of the signal. In addition, by comparing the tSZ 
power spectrum from the resolved clusters with the marginalized 
tSZ power spectrum presented in Sect. 7 we deduce that the mea- 
sured tSZ spectrum includes an additional tSZ contribution from 
unresolved clusters and diffuse hot gas. 




y X 10" 



Fig. 10. ID PDF of the Planck y-map (black) and of the null map 
(blue) for the MILCA method. 



6. Analysis of High Order Statistics 

The power spectrum analysis presented above only provides in- 
formation on the 2-point statistics of the Compton parameter dis- 
tribution over the sky. A full characterization of the field can be 
performed by studying the higher order statistics in the ID PDF 
of the map, or by measuring the bispectrum. 



6.1. 1 D PDF analysis 

We performed an analysis of the ID PDF of the NILC and 
MILCA reconstructed Compton parameter maps. For the tSZ ef- 
fect we expect an asymmetric distribution with a significantly 
positive tail (Rubino-Martln & Sunyaev 2003). We thus focus 
on the asymmetry of the distribution and its unnormalized skew- 



ness. First, we filter the maps in order to enhance the tSZ signal 
with respect to foreground contamination and noise. To avoid 
residual point source ringing effects near the edges of the com- 
bined PSMASK and GALMASK50 masks we apodize them. We 
follow the approach of Wilson et al. (2012) and use a filter in 
harmonic space, constructed from the ratio between the angular 
power spectrum of the expected tSZ signal in the FFP 6 simula- 
tions and the power spectrum of the null y maps. We smooth this 
ratio using a 21 -point square kernel and normalize it to one by 
dividing by its maximum value. Notice that this filter only selects 
the multipole range for which the tSZ signal is large with respect 
to the noise, and thus, it does not modify the non-Gaussianity 
properties. Furthermore, we have found that the filter used here 
behaves better than the more traditionally-used Wiener filter, as 
it is less affected by point-source ringing. Then, the ID PDF of 
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the filtered Compton parameter map, P(y), is computed from the 
histogram of the pixels. 

Figure 9 shows the ID PDF for the FFP 6 simulation maps 
combined using the weights of the MI LCA linear combination 
of the real data. We present in red the ID PDF of the tSZ effect, 
which is clearly asymmetric, with a positive tail as expected. 
Moreover, the asymptotic slope of this red curve at high values 
of y scales almost as P(y) oc y~ 2 - 5 , implying that the underlying 
source counts should scale in the same way (i.e., dn/dY oc Y~ 2 5 ). 
This is the predicted scaling behavior for clusters (e.g., de Luca 
et al. 1995; Rubino-Martm & Sunyaev 2003), and indeed, it is 
the scaling that we find in the actual number counts of clus- 
ters in the simulation used. Similarly, the ID PDF for radio 
sources (black) is also asymmetric, but with a negative tail. By 
contrast, the clustered CIB (blue) and diffuse Galactic emission 
(cyan) distributions are to first approximation symmetric. From 
this analysis we see that the filtering enhances, as expected, the 
tSZ effect with respect to foregrounds and therefore helps their 
discrimination. 

For illustration, Fig. 10 shows the ID PDF for the MILCA 
Compton parameter map in black. This is the convolution 
of the ID PDF of the different components in the map: tSZ 
effect; foregrounds; and noise. Indeed, it clearly shows three 
distinctive contributions: a Gaussian central part that exceeds 
slightly the contribution from noise, as expected from the null 
map ID PDF (cyan curve); a small negative tail, corresponding 
most likely to residual radio sources; and a positive tail cor- 
responding mainly to the tSZ signal. A direct computation of 
the slope of the full P(y) function in Fig. 10 shows that it con- 
verges to -2.5 fory > 10 5 , as predicted from the cluster counts. 

A simple analysis of the measured ID PDF can be done by 
considering the asymmetry of the distribution: 

A= P(y)dy- P(y)dy, (13) 

where y p is the peak value of the normalized distribution 
(J P(y)dy = 1). In addition, the non-Gaussianity of the positive 
tail can be quantified by 

p+oo 

A= [P(y)-G(y)]dy, (14) 

with G(y) the expected distribution if fluctuations were only 
due to noise. For the NILC Compton parameter map we find 
A = 0.185 and A = 0.065. Equivalently, for the MILCA 
Compton parameter map we find A = 0.26 and A = 0.1 1. These 
results are consistent with a positive tail in the ID PDF as 
expected for the tSZ effect. The differences between the NILC 
and MILCA results come mainly from the difference in filtering. 

Alternatively, we can also compute the skewness of the 
obtained distribution, j y 3 P(y)dy/ (J y 2 P(y)dyj . Following 
Wilson et al. (2012) we have chosen here a hybrid approach, by 
computing the unnormalized skewness of the filtered Compton 
parameter maps outside the 50% sky mask. In particular we have 
computed the skewness of the Planck data Compton parameter 
maps /y 3 \, of the null maps (^null)- ^ or tne FFP ^ simulations, 
we computed these for the tSZ component (y FFP6 sz ) and for the 

sum of all astrophysical components (y FFP6 all)- Table 2 shows 
the results for the NILC and MILCA maps. The different filter- 
ing function derived for the NILC and MILCA y-maps prevents 



Table 2. Unnormalized skewness multiplied by 10 18 . 



Method (y 3 ) (yl ULL ) (f FFP6Jiz ) {y 3 FFP6ALL } 

NILC . . . 1.78 -0.0001 2.17 2.09 
MILCA . . 1.50 0.0004 1.46 1.21 



a direct one-to-one comparison of the skewness. However the 
comparison for each map with the FFP6 simulations for the tSZ 
component and for the sum of all components clearly shows the 
minor contribution of the foregrounds in both maps. This allow 
us to argue that the measured skewness value is mainly domi- 
nated by the tSZ signal, as one would expect from Figs. 9 and 
10. By comparing the measured and model skewness we present 
in the Sect. 7.2 illustrative constraints on <x g . 

6.2. Bispectrum 

Since the SZ signal is non-Gaussian, significant statistical in- 
formation is contained in the bispectrum, complementary to the 
power spectrum (Rubino-Martm & Sunyaev 2003; Bhattacharya 
et al. 2012). We therefore compute the bispectrum of the 
NILC and MILCA reconstructed Compton parameter maps. 
The results presented here use the binned bispectrum estimator 
described in Bucher et al. (2010) and Lacasa et al. (2012), 
which is also used for the Planck primordial non-Gaussianity 
analysis (Planck Collaboration XXIV 2013). We mask the maps 
with the combined PSMASK and GALMASK50, remove the 
best-fit monopole and dipole outside the mask, and degrade 
the resolution to N^de = 1024 to reduce computing time. We 
use a multipole bin size A( = 64 and a maximum multipole 
f m „ = 2048 for the analysis. To correct for the bias introduced 
by masking we have produced non-Gaussian simulations with 
a tSZ-like bispectrum and we have convolved the simulated 
maps with a Gaussian beam of 10' FWHM. We compute the 
bispectrum of the simulated full-sky and masked maps and 
measure the average ratio between the two. This ratio is used 
to correct the measured bispectra and flag unreliable (^1,^2,^3) 
configurations for which mask effects are too large to be 
corrected for. 

We checked that foreground residuals do not significantly 
affect the recovered tSZ bispectrum using the FFP 6 simulations 
described previously. In the case of the MILCA reconstructed 
map (more affected by foregrounds), for example, Fig. 11 
shows the tSZ bispectrum as well as the (absolute value of the) 
bispectra of the different foreground residuals, in some special 
configurations, namely equilateral ((, (, (), orthogonal isosceles 
(£,(, y/lf), flat isosceles {(,£,!€) and squeezed ({ min J,t). The 
foreground residuals yield negligible bispectra, at least one 
order of magnitude smaller than the tSZ bispectrum over the 
multipoles of interest. 

In Fig. 12 we compare the tSZ bispectrum measured on 
Planck data, with the tSZ bispectrum of the FFP 6 simulation 
and with the bispectrum of the maps of detected clusters in the 
Planck catalogue presented above. Clusters from the Planck 
catalogue contribute to an important fraction of the measured 
bispectrum, at least 30% on large angular scales and more on 
smaller angular scales; the bispectrum hence also probes the 
unresolved tSZ signal, as was the case for the power spectrum. 
On large angular scales this may be the signature of the clus- 
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Equilateral 



Orthogonal 




Fig. 11. Binned bispectra of the FFP 6 tSZ map, and foreground 
residuals for the MILCA component separation. The black line 
represents the tSZ bispectrum and the red line the clustered CIB. 
In addition, we plot the bispectrum for the Galactic diffuse free- 
free (orange), and the thermal dust (dark blue) emission. 



Table 3. Marginalized bandpowers of the angular power spec- 
trum of the Planck tSZ Compton parameter map (unitless), sta- 
tistical and foreground errors, and best-fit tSZ power spectrum 
and number counts models (unitless). 
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Fig. 12. tSZ measured bispectra for the MILCA (green) and 
NILC (orange) Compton parameter maps compared with the 
FFP 6 tSZ bispectrum (in red) and with the projected map of 
the catalogue of detected clusters (in blue). Uncertainties on the 
measured tSZ power bispectra are given by the dotted lines. 



tering of less massive dark matter halos inside the large-scale 
structure. Alternatively large angular scales may be affected by 
foreground residuals. 



7. Cosmological Interpretation 

7. 1. Power spectrum analysis 

As a measure of structure growth, the tSZ power spectrum can 
provide independent constraints on cosmological parameters and 
potentially improve their precision. As shown by Komatsu & 
Seljak (2002), the power spectrum of the tSZ effect is highly sen- 
sitive to the normalization of the matter power spectrum, com- 
monly parameterized by the rms of the z = mass distribution 
on 8 h~ l Mpc scales, <xg (Komatsu & Seljak 2002), and to the 
total amount of matter Q m . 

Cosmological constraints are obtained from a fit of the NILC 
F/L cross-power spectrum, for the 50% mask, assuming a three- 
component model: tSZ; clustered CIB; and radio and infrared 
point sources. For I > 60, we can reasonably neglect the Galactic 
dust contamination. For £ > 1411 the total signal in the tSZ 
map is dominated by noise. We thus restrict our analysis to the 
multipole range 60 < ( < 1411. The measured power spectrum, 
Cf, is modelled as follows: 

Cf = Cf z (Q m , tr 8 ) + Acib Cf m + A PS (Cf + Cf ad ). (15) 

Here Cf z (Q m , cr 8 ) is the tSZ power spectrum, Cf lB is the clus- 
tered CIB power spectrum and C' R and Cf' dd are the infrared and 
radio source power spectra, respectively. 

Following Sect. 8, the tSZ spectrum is computed using the 
halo model, the Tinker et al. (2008) mass function and the 
Arnaud et al. (2010) universal pressure profile. In particular, 
we use the numerical implementation presented in Taburet et al. 
(2009, 2010, 2011) and integrate in redshift from to 3 and in 
mass from 10 13 to 5 x 10 15 M . Our model allows us to compute 
the tSZ power spectrum at the largest angular scales. It is consis- 
tent with the tSZ spectrum presented in Efstathiou & Migliaccio 
(2012) and used as a template in the CMB cosmological analy- 
sis in Planck Collaboration XV (2013) and Planck Collaboration 
XVI (2013). 
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Fig. 15. Marginalized bandpowers of the Planck tSZ power spectrum with total (statistical plus foreground) uncertainties (red 
points). The red solid line represents the best-fit tSZ power spectrum model. We also show as a blue solid line the best-ft tSZ power- 
spectrum obtained from the analysis of cluster number counts (Planck Collaboration XX 2013). The tSZ power spectrum template 
used in the CMB cosmological analysis (Planck Collaboration XV 2013; Planck Collaboration XVI 2013) is presented as a green 
solid line. 




<r 8 (fi m /0.28) 3/8 ^cm ^ps 

Fig. 13. 2D and ID likelihood distributions for the combination 
of cosmological parameters crg(O m /0.28) 3 ' 2 ^ 81 , and for the fore- 
ground parameters Acm and Aps. We show the 68% and 95.4% 
C.L. contours here. 



Foreground contamination is modelled following Sect. 5.2.2. 
As discussed there the main uncertainties in the residual power- 
spectrum translate into up to 50% uncertainty in the clustered 
CIB and point source amplitudes. We thus allow for a variation 
of the normalization amplitudes for the clustered CIB, Aqb and 
for the point sources, Aps, with Gaussian priors centred on one 
with standard deviation 0.5. 

We have not considered explicitly the expected correlation 
between the tSZ effect and the CIB (Addison et al. 2012). 
However, using the formalism in Addison et al. (2012), we have 
performed simulations of the expected effect and we found that 
to a reasonable level of approximation the shape of the tSZ and 
clustered CIB cross-power spectrum is very similar to that of 
the clustered CIB power spectrum. Therefore, in our simplified 
modeling, the clustered CIB normalization factor, Aqb, also ac- 
counts for this component. 

We assume a Gaussian approximation for the likelihood 
function. Best-fit values and uncertainties are obtained using an 
adapted version of the Cosmo-MC algorithm (Lewis & Bridle 
2002). Only a% and Q m are allowed to vary here. All other cos- 
mological parameters are fixed to their best-fit values as obtained 
in Table 2 of Planck Collaboration XVI (2013). The normaliza- 
tion amplitudes Aqb and A P $, considered as nuisance parame- 
ters, are allowed to vary between and 3. For the range of mul- 
tipoles considered here, the tSZ angular power spectrum varies 
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Fig. 14. Likelihood samples derived from the MCMC chains. 
The points represent pairs of values of O m and cr^. The clus- 
tered CIB amplitude is colour-coded according to the value of 
Acib, from low (blue) to high (red). The black solid line shows 
the theoretical degeneracy between the two cosmological param- 
eters. 



like C( oc crjj'Q^ 2 . The results are thus presented in terms of this 
parameter combination. 

Figure 13 presents the 2D and ID likelihood distributions 
for the cosmological parameter combination cr 8 Q^ 2 ' 81 , and for 
the foreground nuisance parameters. The best-fit values and er- 
ror bars for each parameter are given by cr 8 (Q m /0.28) 3 2 ^ 81 = 
0.784 ± 0.016, cr 8 = 0.74 ± 0.06, fi m = 0.33 ± 0.06, A aB = 
0.55 ± 0.26 and A PS = 0.14 ± 0.13. It is worth noting that these 
values are obtained in a specific framework: all other cosmolog- 
ical parameters fixed and a fiducial fixed model for the signals. 
Opening up the possibility for relaxing this framework is likely 
to weaken the constraints. 

Figure 14 shows the degeneracy between the two cosmo- 
logical parameters from the MCMC chains, as well as the the- 
oretical degeneracy (solid black line). It also shows the de- 
pendency on Acm (colour coded from low values in blue to 
high values in red). While the combination o"8(f2 m /0.28) 3 2 ^ 81 
is well determined, marginalized constraints on <x 8 and Q m are 
weaker. To check the robustness of our results, we performed the 
same cosmological analysis using the NILC-MILCA F/L cross- 
power spectrum presented in Fig. 6. Although the foreground 
level is different, we find compatible results at the 1 cr level. 
Furthermore, our constraints are in good agreement with those 
derived from the Planck cluster number count analysis Planck 
Collaboration XX (2013), which shows similar <x 8 -Q m degen- 
eracy line. Conversely, our findings exhibit some tension with 
the constraints derived from the Planck primary CMB analysis 
(Planck Collaboration XVI 2013), which finds larger values of 
erg and D. m . As discussed in Planck Collaboration XX (2013), 
the constraints from the SZ signal depend significantly on the 
assumed value of the mass bias. 

The red points in Fig. 15 correspond to the marginalized 
Planck tSZ power spectrum (from the NILC F/L cross-power 
spectrum), compared to the best-fit theoretical model presented 
above (solid red line). Foreground uncertainties are derived from 




e 



Fig. 16. Comparison of the best-fit Planck tSZ power spec- 
trum model with existing models in the literature. We show the 
Planck tSZ power spectrum best-fit models derived in this pa- 
per (red) and in the analysis of cluster number counts (Planck 
Collaboration XX 2013, blue). The tSZ power spectrum 
template used in the CMB cosmological analysis (Planck 
Collaboration XV 2013; Planck Collaboration XVI 2013, green) 
is also shown. We also show the tSZ power spectrum models 
from hydrodynamics simulations (Battaglia et al. 2012, brown), 
from N-body simulations plus semi-analytical dust gas models 
(Trac et al. 201 1, cyan and purple), and from analytical calcula- 
tions (Shaw et al. 2010, black; Komatsu & Seljak 2002, orange). 



the likelihood curves of the nuisance parameters and added in 
quadrature to the statistical uncertainties, providing the total er- 
rors plotted here. Table 3 presents the Planck marginalized tSZ 
power spectrum, statistical and foreground uncertainties, and the 
best-fit tSZ power spectrum model. In the range ^=60-141 1, the 
Planck tSZ power spectrum can be approximated by a power law 
of the form 

{{€ + l)Ce/2n = (1.0 ± 0.2) ^°- 912±0 - 03I > 10 15 . 

The measured tSZ power-spectrum is in remarkable agree- 
ment with the tSZ power-spectrum (blue solid line) computed 
using cluster count best-fit parameters (Planck Collaboration XX 
2013). We also show in Fig. 15 (green line) the tSZ template used 
in the Planck CMB analysis Planck Collaboration XVI (2013). 
This template is renormalized to the measured Planck tSZ power 
spectrum at £ — 959. The difference in shapes of the two spectra 
is due to the different assumptions on the scaling relation be- 
tween SZ signal and mass (see, Efstathiou & Migliaccio 2012). 
We also show the SPT (Reichardt et al. 2012, orange point) and 
ACT (Sievers et al. 2013, cyan point) constraints on the tSZ 
power spectrum at £ — 3000, which are consistent with our best- 
fit model within +2cr and illustrate that the tSZ spectrum starts 
to turn over at higher I. 

In Fig. 16, we compare the Planck tSZ measurements of the 
power spectrum to a set of predicted spectra. We consider the 
models from hydrodynamical simulations (Battaglia et al. 2010; 
Battaglia et al. 2012, brown), from N-body simulations plus 
semi-analytical models (Trac et al. 2011, cyan and purple) and 
from analytical calculations (Shaw et al. 2010, black); (Komatsu 
& Seljak 2002, orange). These models were computed originally 
for the set of cosmological parameters in Hinshaw et al. (2012) 
with cr 8 = 0.8 and have been rescaled to our best-fit value for 
cr 8 'Q^ 2 . We first note the large dispersion in the predicted am- 
plitudes and shapes of tSZ spectrum. These differences reflect 
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the range of methodologies and assumptions used. This of course 
includes the redshift ranges and most of all the mass intervals 
probed by the limited-sizes of the simulation boxes of hydro- 
dynamical simulations. On the other hand analytical predictions 
are sensitive to the model ingredients (mass function, mass bias 
and scaling relations). The Planck data allow us for the first time 
to probe the whole range of angular scales for the tSZ spectrum; 
furthermore and as shown in Planck Collaboration XXIX (2013), 
Planck is particularly sensitive to the SZ signal from massive 
clusters not probed by other experiments. The Planck tSZ mea- 
surement will hence permit in the future to build a coherent un- 
derstanding from the integrated tSZ contribution of the whole 
population of clusters, including resolved and nearby clusters, 
the correlated SZ signal and possible diffuse hot gas. 

7.2. High order statistics 

The estimates of tSZ non-Gaussianity, e.g., unnormalized skew- 
ness and bispectrum, are very sensitive to <x 8 . Using the mod- 
els presented in Sect. 3 we can show that the unnormalized 
skewness of the tSZ fluctuation, ^r 3 (n)^ scales approximately 

as o"g , whereas the amplitude of the bispectrum scales as 
(j? with a = 1 1-12, as shown by Bhattacharya et al. (2012). 

We derive constraints on <x 8 by comparing the measured 
unnormalized skewness and bispectrum amplitudes to those 
obtained from simulations of the tSZ effect. This approach 
is strongly limited by systematic uncertainties and the details 
of the theoretical modelling (see Hill & Sherwin 2013, for a 
review). 

From the measured unnormalized skewness of the filtered 
MILCA and NILC Compton parameter maps discussed in 
Sect. 6. 1 and by comparing them to the value measured in the 
FFP 6 simulations we can derive constraints on cr 8 . Uncertainties 
due to foreground contamination are computed using the FFP 6 
simulations and are accounted for in the final error bars. The 
tSZ component of the FFP 6 simulations was obtained from 
a hybrid simulation including a hydrodynamic component for 
Z < 0.3 plus extra individual clusters at z > 0.3, and with 
cr 8 = 0.789. Using these simulations we obtain <x 8 = 0.775 for 
NILC and o- 8 = 0.783 for MILCA. Combining the two results 
and considering model and foreground uncertainties we obtain 
cr 8 = 0.779 ± 0.015 (68 % C.L.). Notice that the uncertainties 
are mainly dominated by foreground contamination. Model 
uncertainties here only account for the expected dependence 
of the unnormalized skewness with <x 8 , as shown in Sect. 3. 
We have neglected here, as was also the case in Wilson et al. 
(2012), the dependence on other cosmological parameters. We 
have also not considered uncertainties on the combination of the 
hydrodynamical and individual cluster simulations. Because of 
these constraints, our error bars might be underestimated. 

The comparison of the measured bispectrum on the Planck 
Compton parameter maps with the FFP 6 simulation tSZ bispec- 
trum shows an offset of about a factor of two on small angular 
scales, 300 < £ < 700, that we attribute to the differences in 
cosmological parameters. Using the scaling of the bispectrum 
with <r 8 , its uncertainty, as well as the uncertainty on the bis- 
pectra ratio, we obtain cr 8 = 0.74 ± 0.04 (68 % C.L.). As for 
the unnormalized skewness, we neglected here the dependence 
on other cosmological parameters and the uncertainties on the 
FFP 6 simulations, and thus, the error bar might be again some- 
what underestimated. 



8. Conclusion 

Because of its wide frequency coverage from 30 to 857 GHz, the 
Planck satellite mission is particularly well suited for the mea- 
surement of the thermal Sunyaev-Zeldovich effect. Taking ad- 
vantage of this, and working with the Planck frequency channel 
maps from 100 to 857 GHz, we have reconstructed the tSZ signal 
over the full sky with tailored component separation methods. In 
this paper, we presented the first all-sky tSZ map quantified in 
terms of the Compton parameter and with an angular resolution 
of 10'. 

We have characterized the reconstructed Planck all-sky 
Compton parameter map in terms of blind detection of tSZ 
sources, and the angular power spectrum and higher order statis- 
tics via the study of its ID PDF and bispectrum. In all cases we 
have identified, characterized and carefully modelled the con- 
tamination by foreground emission. These are mainly due to 
diffuse Galactic thermal dust emission at large angular scales 
(t < 60), and clustered CIB and Poisson-distributed radio 
and infrared sources at smaller angular scales (dominating at 
i > 500). Diffuse Galactic thermal dust emission is tackled 
via a conservative masking of the brightest 50% of the sky in 
the Planck 857 GHz channel map. The CIB and point source 
contamination are modelled in agreement with the findings of 
Planck Collaboration XVIII (2011) and Planck Collaboration 
(2013). 

We have produced the first measurement of the SZ power 
spectrum on large angular scales ranging from 0.17° < 9 < 3.0°. 
In this range, the tSZ power spectrum is almost insensitive to the 
physics of cluster cores. The detected tSZ signal likely arises 
from the contribution of warm and hot diffuse gas distributed 
within groups and clusters, sampling the whole halo mass func- 
tion, as well as within the larger scale filamentary structures. 

We have modelled the tSZ power spectrum via a halo-model 
analytical approach, in order investigate its dependence in terms 
of <x 8 and Q. m and to test it against the measured Planck tSZ 
power spectrum. Moreover, we performed an analysis of the ID 
PDF and bispectrum of the Planck y-map to infer independent 
constraints. We find, in the present framework, that the best-fit 
normalization parameter cr 8 from the three independent analysis 
ranges between (0.74 + 0.06) and (0.78 + 0.02) at 68% C.L. 

These values are lower than the ones derived from analysis 
of primary CMB anisotropics (Planck Collaboration XVI 2013). 
More refined analysis and modelling will be needed to under- 
stand this difference. This tension may have different possible 
origins. Some of the difference may be due to the hypotheses in 
the tSZ modelling, e.g., mass bias (see Planck Collaboration XX 
(2013) for a detailed discussion on its effect on cluster counts). 
Other differences could arise from the foreground modelling in 
particular at high frequencies, above 217 GHz. 

The observed consistency between constraints derived from 
the cluster number counts in Planck Collaboration XX (2013) 
and from the present work seems to provide a coherent view of 
the gas content in halos and in larger scale structures. As such, 
this Planck tSZ measurement constitutes the first step towards 
building a coherent understanding of the integrated tSZ effect 
due to cosmic structure on all scales and density contrasts. 
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Appendix A: Foreground contamination in the 
final tSZ power spectrum 

Since we are using modified Internal Linear Combination 
methods to estimate the final Planck Compton parameter map 
we can write it as 

y{6, 4>) = Yj Yi W v {e ' ^ ( Fh(e > & * Mv(e > 0) ) ' (A ' 1} 

v b 

where M v {8, <p) is the Planck map for frequency channel v, 
F h (9, </>) is a circular filtering function for the multipole interval 
b and Wy (9, 0) are the weights of the internal linear 
combination in that multipole range. Decomposing y(8, <p) in 
spherical harmonics we obtain 

%m = J] J] W b Y(M F> * M n , w , . (A.2) 

v b 

Then using spherical harmonic convolution properties (see for 
example Tristram et al. 2005) and assuming overlap in the 
multipole range selected by the filter functions, F b r then the 
power spectrum is given by 

Cf = ZZZZZ<^cf^. (A.3) 

b b' v V I' 
W b ,W b ', 

where M ( \, v represents the mode coupling matrix associated 
withW^x W%. 

For each Planck channel the sky signal can be expressed as the 
sum of multipole components including CMB, tSZ, diffuse 
Galactic emission, radio and IR point sources and clustered 
CIB, such that the Planck Compton parameter is given by 

y = y + y CMB + Yjy c > (AA) 

c 

where c sums along the different foreground contributions. By 
construction ycivtB = and thus, assuming no correlation 
between foreground components, the estimated tSZ spectrum 
can be expressed as 

d/ = Cf z + Ycf' f . (A.5) 
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Using expression A. 1 we write 

y\6, cf>) = J] X ® ^ * M ^ *>) ' (A - 6) 

v 6 

and thus using equation A. 3 we have 

i i' v v' f 



The latter expression can be simplified assuming a common 
spatial distribution of the foreground emission across 
frequencies and a well defined spectral energy density, f£, so 
that it reads 

C V -XXXX m ct mc - (a.8) 

b v v' (' 

Let us now look at the cross-correlation between the estimated 
Compton parameter map and a particular sky component at one 
of the observation frequencies M V C (Q, <p) - Y*tm T v . ( YttJfi, <p). 
We define the cross-power spectrum as 



(A.9) 



m=-{ 

and the statistical expectation of this quantity reads 

m——l by I'm' 

l 1 

— YjJuYu W ^ F 1' C % 6 «' S ^ (A.ll) 

(A. 12) 



It- 



m=—£ by I'm' 



2C 



1 V V w bV p b r w 

+ j /j/j n U;mm r 1^ c;C ' 



m=-l by 



Assuming we have a measure of C v c v ( and a way to derive 
W%! F%, we can compare the measured cross-correlation y xc 
to its theoretical expectation and thus have a consistency check 
on the component contribution to y. 
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